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from a gap equation which is either obtained from the Nambu Jona-Lasinio 
' Lagrangian or from a fit to observables calculated on the lattice. The latter 
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\ shows an effective confinement. At thermal equilibrium the stability is stud- 
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' ied within thermodynamics (mechanical stability) and via a linear response 

analysis of the Vlasov equation. The instabilities related to a first-order phase 
qq | transition are found. For a plasma expanding in three and one dimensions far 

from equilibrium a new type of instability, called dynamical, appears. The 
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| relation to cluster formation is shown in a molecular dynamics calculation. 
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I. INTRODUCTION 



The experimental search for the quark-gluon plasma in ultra-relativistic nuclear collisions 
requires the development of reliable and predictive models of the formation of the plasma 
and of its hadronization in order to be able to extract observables which signal the existence 
of a parton plasma in the high energy density region of the colliding nuclei. While lattice 
calculations have yielded important information on the equilibrium properties of a parton 
plasma (at least for baryon chemical potential /z = 0), they are unable yet to describe the 
space-time development of the plasma formation and its subsequent decay into hadrons. The 
formation phase of the plasma is calculated usually in cascade models [Q, which predict high 
values for the energy densities in the center of the collision. These values should be sufficient, 
according to lattice QCD calculations, for a plasma of deconfined quarks and gluons to be 
formed. A satisfactory microscopic calculation of the hadronization process in space-time 
is not yet available, in particular because of the complexity of the confinement problem. In 
this situation model studies are appropriate to elucidate the physics of the process. The 
present work is one of them. 

In our work the partonic phase will be described within a quasi-particle model with an 
effective mass. The simplicity of the model allows to describe at the same time the equi- 
librium properties and the expansion dynamics. No hadronization is yet included. The 
medium dependent mass m (in the case of thermal equilibrium m(T, //), and m(x,t) for an 
expanding parton gas) is calculated self-consistently from a gap equation which relates the 
mass m to the local parton distribution function in phase space. The self-consistent mass 
itself determines the momentum distribution in the case of thermodynamical equilibrium or 
the evolution equation (the Vlasov equation) of the phase-space distribution in the nonequi- 
librium situation. We study two quasi-particle models, which are in a sense complementary: 

1. The gap equation is the one of the Nambu Jona-Lasinio (NJL) model |2|,[|] with non- 
zero current quark mass. For /x = it shows a cross-over transition around a "critical" 
temperature T c and a first one at low temperatures and finite baryon density. 
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2. We introduce a phenomenological parton model (PPM): Its gap equation is adjusted 
such that the quasi-particle model reproduces the lattice data for the energy density 
e(T) and the chiral condensate (^/^/){T) around T c (for \i = 0) j4|. It is described in 
Ref. [[[J and generalized to p, ^ 0. It shows effective confinement in that the system 
cannot expand indefinitely, because the effective parton masses would grow to infinity. 

For these two quasi-particle models we study the stability of the system in three physical 
situations. 

• Within equilibrium thermodynamics the region in the plane of temperature T and 
baryon density ps is calculated where dP/dps < 0, i.e., the system is mechanically 
unstable. At the same time the region of the mixed phase in the T-pg plane is calcu- 
lated which indicates the presence of a first order phase transition. 

• Using the method of linear response, the stability of the space-time dependent Vlasov 
equation around thermal equilibrium is investigated and the exponential growth rates 
are determined for the region of instability. 

• The stability analysis of the Vlasov equation is also applied to situations far from equi- 
librium as they occur in three-dimensional (Hubble) and one-dimensional (Bjorken) 
expansions of a system of partons. 

• In order to show that the instabilities are related to hadronization, we have studied the 
expansion and clusterization of an initially compact system of partons via molecular 
dynamics. 

II. THE MODEL 

The basic assumption of our work is the quasi-particle approximation for a gas of "par- 
tons" (quarks and gluons) with an effective medium dependent mass The thermody- 
namic potential of the quasi-particle system composed of fermions, antifermions and gluons 
with degeneracies gf, gj = gj and g g , respectively, is 

3 



p 



d 3 p p 2 



9fff(p) + 9fff(P) + 9gfg(p) 



- V(m) , 



(1) 



(2tt) 3 3E 

where fa is a Fermi or Bose distribution function and V(m) is the mass dependent bag 
pressure. In our calculation we use the same effective mass for quarks and gluons. The 
parton mass depends on the medium through the gap equation : 
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which follows from the thermodynamical relation 
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Eq. ([]]) is not the most general expression for the pressure. For instance the bag pressure 
could also depend on the baryon density. Since lattice data are only available for zero baryon 
density, we cannot explore this possibility. For the NJL model, which has no gluons, the 
explicit form for dV/dm can be derived from the NJL Lagrangian within the mean field 
approximation. One has 

(m — mo) f d 3 p m 



dV 
dm 
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where G is the four fermion coupling constant, A a three-momentum cutoff, and mo the 
current quark mass. At zero baryon density the NJL model predicts a second order phase 
transition (for mo = 0) or a crossover behavior (for mo > 0) (see Fig. |1]). At finite baryon 
density and low temperature the systems shows a first order phase transition for moderate 
values of the current quark mass (Fig. Q). We use the parameters m = 2 MeV, A = 653 MeV 
and GA 2 = 2.1. 

For the second model, called the phenomenological parton model (PPM) [||, the tem- 
perature dependence of the parton mass m(T) is obtained by fitting the energy density 
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to lattice QCD data around T c . The parton mass in the high temperature region T ^> T c 
is required to go smoothly towards zero in order to reproduce qualitatively the vanishing 



chiral condensate at high temperature. For the present study we use the same procedure to 
obtain m(T) as explained in Ref. ||, except that here we use Fermi and Bose distributions 
for the fermions and gluons, respectively, and take gf = 12 corresponding to the two-flavor 
QCD. The energy scale is set by the critical temperature T c = 140 MeV for the two-flavor 
QCD. Fig. [3] shows m(T). Once it is known, one can use Eq. (|2|) to obtain dV/dm and V as 
a function of m. As explained in the beginning of this section, the quasi-particle model is 
completely defined once the bag pressure V{m) is given. Fig. [| compares the bag pressure 
V(m) of the PPM model with the one of the NJL model. Note that we use gf = gf = 12 
in both models for the fermions, but for the gluons g g = in the NJL and g g = 16 in the 
PPM model. The use of the Fermi and Bose distribution functions and the assumption 
that V depends only on m allow one to explore the parton mass and the thermodynamical 
properties of the system in the whole plane of temperature and baryon density, even where 
there are no lattice data. Fig. |5] shows the parton mass as a function of the chemical potential 
at T = 100 MeV for the phenomenological parton model. A first order phase transition 
is present at finite baryon density although at zero baryon density we have, following the 
lattice QCD, only a crossover behavior around T c . This behavior, crossover transition at 
/i = and first order one at finite /i, is also found for the NJL model. 

The quasi-particle models, whose thermodynamical properties are described above, can 
also be used to describe the space-time evolution of parton system. The underlying equation 
is the Vlasov one 

d t fi(t,x,p) + — ^ - ■ V x fi(t,x,p) - E(tx p) Vxm ^^ ' V r>hit^,P) = , (6) 
where the quark, anti-quark or gluon distribution functions fi at time t, at space point 



x and at three momentum p describe on-shell partons with energy E = y/m 2 + p 2 . In a 
quasi-particle model there are no interactions between the quasi-particles, and therefore 
there are no collision terms on the r.h.s. of Eq. (|6]). The space-time dependent mass m(t, x) 
is calculated self-consistently from the gap equation as a function of the non-equilibrium 
parton distribution functions fi(t,x,p) : 



dV ^ r d 3 p m(t,x) 

Eqs. (|6]) and (0) have to be solved simutaneously for an evolving parton gas. This has 
been done already for the phenomenological parton model with the result that during the 
expansion the partons acquire high masses and the expansion stops |J . We call this behavior 
"effective confinement" . 



III. MECHANICAL STABILITY AND FIRST ORDER PHASE TRANSITION 

For a system at thermodynamic equilibrium the mechanical instability is defined by the 
condition 

(£),<■■ 

With the gap equation @ the condition (|8|) is equivalent to 

(£),<"■ 

The negative derivative of the pressure or the chemical potential with respect to the 
density signals an instability of the homogeneous state of the system. Generally such an 
instability is related to a first order phase transition || . 

We discuss the numerical results for the two models, and begin with the NJL model. 
The pressure as a function of the baryon density is shown in Fig. ^ for several values of 
the temperature. For T = for example the pressure is not a monotonic function of the 
density and a first order phase transition occurs. For the interval B C in Fig. |6] the system is 
unstable since one has OP/dps < 0. The thermodynamical system is in a metastable state 
in the regions (intervals A-B and C-D in Fig. ^ around the unstable region. Performing 
a Maxwell construction amounts to choosing the correct solution of the gap equation for 
certain values of fi where the gap equation has 2 or 3 solutions for the mass (short dashed 
line in Fig. ^]). The so constructed mass as a function of the chemical potential has a 
discontinuity at the critical chemical potential \x c (see Fig. The energy density and the 
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baryon density are also discontinuous at \i c . The discontinuity in ps means that the whole 
region AD of the baryon densities is excluded from the phase diagram. This region is called 
mixed phase and is shown in Fig. |7| in the plane of temperature and baryon density. The 
region of mechanical instability defined by Eq. (^) is a part of the mixed phase. Inside the 
region of mechanical instability the system is unstable with respect to non-homogeneous 
perturbations. In the whole region of mixed phase it is preferable for the system to separate 
into zones belonging to the low or high density thermodynamically stable phases. 

For the phenomenological parton model, Fig. || shows the pressure as a function of the 
baryon density for several temperatures. The behavior is qualitatively similar to the NJL 
model: For T < 120 MeV the pressure is not a monotonic function of the density, and a 
Maxwell construction is required. However, quantitatively the situation for the PPM model 
is rather different and we have chosen the insert to show the situation for small p#. 

The regions of mechanical instability and metastability in the T-ps plane are shown in 
Fig. ^| for the phenomenological parton model. A region with temperatures below 120 MeV 
and finite baryon density is unstable. As expected, for large temperature and/or baryon 
density the system is stable against mechanical perturbation, which is natural since then 
the partons are almost massless. However, the region around the y-axis (pg = and T > 0) 
is also stable. 

The phenomenological parton model becomes questionable at values of T and pB where 
the parton masses become much larger than the hadron masses, since hadronization is not 
accounted for. This is the case for the interior of the region of the mixed phase. However, 
close to the boundary of this region the model should correctly describe the thermodynamic 
situation including the effect of confinement. 

Let us imagine a heavy ion collision and suppose that the system evolves at or close to 
thermal equilibrium: A region of deconfined parton plasma is created at high temperature 
and high baryon density in the stable regions of the phase diagrams in Figs. ^ or |^. During 
its expansion, the system reaches the region in the T-ps plane where the local homogeneous 
parton distribution is unstable. In the spinodally unstable region any perturbation of the 
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homogeneous distribution amplifies and the system breaks up into regions ( "clusters" ) with 
larger density. This is the onset of hadronization. However, our quasi-particle models are 
not realistic in these regions, their and the structure, sizes and masses of the clusters cannot 
be reliably calculated. Sect. |V11| presents an illustration of the clusterization process. 



IV. LINEAR INSTABILITY OF THE VLASOV EQUATION NEAR THERMAL 

EQUILIBRIUM 



In the semi-classical limit the collisionless expansion of the plasma can be described by 
the Vlasov equations @ together with the gap equation (|7|). In this section we address the 
linear stability problem of the Vlasov equation for the quasi-particle parton model around 
the equilibrium solution. We consider a small perturbation 8 fa around the homogeneous 
solution fa of the Vlasov equation : 



fi(t,x,p) = fi(p) + 5fi(t,x,p) 



(10) 



where the index i denotes the fermion, anti-fermion or gluon distribution. As it is standard 
practice f^Q, one linearizes the Vlasov equation and takes the Fourier transform in the t 
and x variables to obtain 
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where v = p/E, and 5m(k,uj) is the change of the mass induced by the perturbation of the 
distribution functions 5 fa, and has to satisfy the linearized gap equation 
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The fluctuation of the phase-space distributions 5fa(ui,k,p) can be expressed by the mass 
fluctuation 5m(k, u) using the linearized Vlasov equation (|TT| ) 
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and is inserted into the linearized gap equation: 
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The propagation of the disturbance 5m is only possible, 5m ^ 0, if u and are related, 
u(k), such that the term in the square brackets is zero : 
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The dispersion relation uj(k) defined by Eq. ( |T6"D has two general properties. The solutions 
for the frequency u come in pairs ±u. The frequency is proportional to the wave-vector 
Uk ~ k, which is a general property as long as there is no scale in the dynamical equations. 
In our models the absence of scale is due to the fact that the relation between the density 
and the mass is local (zero-range interaction). 

In the case of an imaginary solution of Eq. (|16D for the frequency, u>k = iT^, the initial 
disturbance grows exponentially in time 



5f(x,t,p) ~ exp(F k t) 



(18) 



i.e. an instability appears. As a consequence the system breaks up into smaller pieces. The 
fact that the instability rate is proportional to the wave-vector (T^ oc k) means that the 
instability rates are unbounded, becoming arbitrary large for large k or small sizes. On the 
other hand we expect that the true clusters into which a low density parton system breaks 
up are the hadrons. The preference of the nature for the two or three quark clusters cannot 
be described in our simple model without color and with only local mean field interaction. 
A finite range of interaction or a minimal size of the cluster would regularize this unphysical 
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behavior. In that case there exists a maximal growth rate of instabilities for a specific 
wave- vector k [[/J . We note that the presence of higher gradients in the dynamical equations 
would change the dispersion relation and limit the growth of excitations with large k ||. 
In our case, the wave-vector for the physically most interesting instability, namely hadron 
formation, should be of the order of 1.0 fm" 1 . 

In the case of the stability analysis around thermal equilibrium, the instability observed 
in the Vlasov equation coincides with those of the mechanical instability and is thus related 



to the occurrence of a first order phase transition. Fig. [10| shows the contours of the growth 
rate F/k in the T-ps plane for the phenomenological parton model. The region of F/k > 
corresponds to the mechanical instability region studied in the previous section. The analysis 
based on the Vlasov equation thus reproduces the results of the thermodynamic analysis but 
goes further in that it also gives values of the growth rates in the unstable region. In the 
instability region close to the phase boundary, where we trust the model, calculated growth 
rates are of the order F/k = 0.1. For k = 1.0 fm _1 the instability growth time 1/r is of 
order 10 fm/c and rather small on the time scale of nuclear reactions. 



V. A GENERAL RESULT FOR THE CASE OF ZERO BARYON DENSITY 



In this subsection we analyze the situation of finite T and pb = 0. We prove that for any 
gap equation which yields dm/dT < 0, the Vlasov equation is stable against perturbations 
of the homogeneous equilibrium distribution. For T > and pb = the dependence m(T) 
of the mass on the temperature is obtained from the equilibrium gap equation 

dV x , f d 3 p m 
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where Fi(E/T) is the equilibrium distributions for p = written as a function of the variable 
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Inserting Eq. (pO|) into the relation (|T6|) a dispersion relation for u and A; is obtained which 
depends on the derivative dm/dT : 




For any thermal distribution F(E/T) its derivative F' (E/T) is always nonzero and negative. 
The solutions for the collective frequency u depend on the sign of dm/dT. Using the 
properties of the function xi^/kv) for real and imaginary uj || one obtains 

• In the case of dm/dT > 0, an imaginary solution of the dispersion cannot be excluded. 
However, no solution has been found for a version of the phenomenological parton 
model studied in Ref. pfl, where the thermal mass increases at high temperatures. 
The absence of a solution can be traced to the fact that the thermal mass increases at 
high temperature only like 

dm m , . 

If < r < 22 > 

• In the case of dm/dT < 0, there is no real or imaginary solution of the dispersion 
relation, and therefore the system is stable. This case applies to the two versions of 
the quasi-particle models studied in this paper, since dm/dT < (Figs. |l]and|3|). Any 
initial perturbation of the homogeneous distribution will be Landau damped and any 
oscillations or instabilities do not propagate or grow. 

Note that the discussion is quite general. It applies to any effective quasiparticle theory 
of the type described in Sect. f|, independently on the statistics of the particles and on the 
specific dependence of the mass on the temperature. In particular we confirm the results 
of Ref. |H| for the Nambu-Jona-Lasinio model at finite temperature and ps = 0, where no 
instabilities have been found. 

VI. LINEAR INSTABILITIES IN AN EXPANDING SYSTEM 



In the previous sections the instabilities of the Vlasov equation are studied at thermal 
equilibrium within the NJL and the phenomenological parton models. While this analysis is 
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important in order to have a bridge to the condition of mechanical instability in equilibrium 
thermodynamics, it may not apply to the case of an expanding plasma, since situations 
far from equilibrium may occur, if thermalization is too slow. In the following we neglect 
thermalization altogether. We study the instabilities of a collisionless expanding plasma 
described by the Vlasov equation. The expansion of a thermal system quickly develops 
into a non-thermal distribution for which the instabilities may be of different nature. They 
appear also for ps = 0. We study two examples of solutions for the Vlasov equation, the 
Hubble and Bjorken scenarios. 



A. Three-dimensional expansion 

The solution for an expansion in three dimensions, which is studied in this subsection, 
is called the Hubble solution, since it describes a homogeneous expansion where the local 
expansion velocity is proportional to the space-time vector x^. For the three-dimensional 
expansion of the plasma any function 

f(p,x,t) = f (s) (23) 

of the variable 

8 = -J(p^y-m(ry (24) 

TO V 

solves the Vlasov equation, and the mass depends only on the proper time r = \/t 2 — x 2 
and Ufj. = x^/t is the local expansion velocity [|n[]. If the system is in thermal equilibrium 
at the initial proper time r with mass m(ro) then the solution can be written in the form 

f{p, x, t) = F^s 2 + m(T y/T) , (25) 

where F(x) = exp(— x) is the thermal distribution function (In this section we restrict 
ourselves only to Boltzmann statistics, since we consider the case /i = 0. We have checked 
that the stability condition and the instability growth rates are not much affected by the 
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statistics.). For the case m(ro) = 0, / is only a function of s/T for all space-time points. In 
the local rest frame x = 0, r = t one has s = \p\r/r and 

f(p,r) = F(\p\r/Tr ) (26) 

for all times. This distribution is a thermal one with time-dependent temperature T{r) = 
Ttq/t, provided the mass remains zero. However the self-consistent mass m(r) calculated 
from the gap equation 

increases rapidly with time r (Fig. [0]). In this case F(pr /Ttq) is no more a thermal 
distribution. In order to analyze the stability of the solution (p6| ) we linearize the Vlasov 
and the gap ( p7|) equations. The procedure leads to a complicated expression, since both 
the zero order solution F(s/T) and the linear perturbation 5f depend on space and time. 
If we restrict ourselves to the vicinity of the point x — 0, the linearized Vlasov equation 
simplifies, since the x-dependence of the zero order solution is weak: 

d t 5f(p, x ,t) + -^-- VJf{p, x, t) - ^}flvjm(x, t) ■ V p F(pr/Tr ) = , (28) 
E t) E[t) 



where now r = t. If one wants to proceed as in Sect. [TV] and introduces the Fourier transforms 
of Sf and 8m with respect to x and t, one has to neglect the time dependence of the zero 
order solution F(pr/Tr ) and m{r) in comparison to the expected stronger time variation 
of the perturbation 5f(p, x, t). This is equivalent to an adiabatic approximation considering 
t as an external parameter which changes the properties of the medium. After Fourier 
transformation one obtains the familiar form of the linearized Vlasov equation 

(u — k- v)5f(p, k, uj) + —k ■ vd p F( P T/Tr )5m = . (29) 



Together with the gap equation fl2?| ) one gets the dispersion relation in the form: 

/;|^V/r4(£)- P ^)=o, (so, 
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which has an imaginary solution for the frequency if 

r>-'™( 1 - P ^)>o. (31 > 

The dispersion relation (^) determines the ratio u/k. For > one has an exponen- 
tially growing instability. Fig. [12] shows the growth rate of the instability as a function of 
the expansion time r/r . If we chose a typical value k = 1.0 fin -1 , the growth rate refers 
to a spinodal instability of the size 1 fm. At t/tq = 1.9 the system becomes unstable and 
quickly reaches a regime of rather short instability time scales of order 1/T = 1.0 fm/c. 

A necessary condition for the appearance of the linear instabilities in the expanding 
plasma is that the mass increases significantly around T c . For example we have found that 
the linear instabilities do not appear for the NJL model or for the three-dimensional Hubble 
expansion in the phenomenological parton model where the parton mass in vacuum is taken 
smaller than ~ 3T C . Note that an instability develops for the case pB = 0, although there is 
no first order phase transition as discussed above. 

B. One-dimensional expansion 

For a one dimensional homogeneous expansion of the plasma (Bjorken scenario) the 
Vlasov equation is solved by any function of the scaling variable s and of the transverse 
momentum p± [[□]] 

f(p,x,t) = f (s, P± ) , (32) 

with 

s = -vV%) 2 -pi-m(r) 2 (33) 

and = (t,0,0, z)/t. In the local rest frame (z = 0) and for a thermal initial state with 
m i T o) — the one-dimensional solution takes the form 

f(p, x, t) = F(J(T Pll /T y+pl/T) = F(E T /T) . (34) 
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The mass m(r) can be calculated from the gap equation and is shown in Fig. [11]. In the 
one-dimensional case the mass increases more slowly than in the three-dimensional case since 
only the longitudinal momentum is dilated. 

Linearizing around the one- dimensional homogeneous expansion, we consider two kinds 
of perturbation, longitudinal Sf(p,z,t) and transverse Sf(p,x±,t) ones. Any general linear 
solution of the Vlasov equation can be factorized into the transverse and the longitudinal 
ones. The approximations in deriving the dispersion relation from the linearized Vlasov and 
gap equations are the same as in the three-dimensional case. One obtains 



P±dpA 



dp\\ 



F\E T /T) 
EE^ 



Vll , p 2 d In t/t 
cu — k\\v\\ " mdm 







(35) 



where v» = p\\/E and k» is the longitudinal wave- vector of the perturbation. The above 
equation has a purely imaginary solution uo = ±iT if 

2 d In t /tq 



r F (E T /T) 
p±dp ± / dp\\ — — 

J~co EjI1i t 



1 - pi 



mdm 



> 



(36) 



In the transverse direction the dispersion relation reads 



oo r2ir 

P±dp± / . 
o Jo 



dp\ 



F\E r /T) 
EE T 



k±v±cos(p ^ 2^( r / T b) 



dm 2 
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uo — k±v± cos 

where k± is the transverse wave- vector of the perturbation. The transverse dispersion rela- 
tion has a purely imaginary frequency solution if 



P±dpA 



dp\\ 



F\E T /T) 
EE T 



1 — p] 



■d(r/r f 



dm 2 



> . 



(38) 



Comparing the conditions for the longitudinal and transverse instabilities (Eqs. fl36|) and 
(|38|)) one finds that the system becomes unstable first in the longitudinal direction, since 

d(r/r ) 2 dlnr/ro 



dm 2 



> 



mdm 



(39) 



Furthermore the longitudinal instability growth rate is larger than the transverse one. Fig. |12 
shows the time dependence of the longitudinal instability growth rate for the one-dimensional 
expansion. The growth rate T is much smaller than in the three-dimensional case and is 
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comparable to the growth rates in equilibrium at finite baryon density. For the wave- vector 
k = 1.0 fm/c one has a maximum value of T ~ 0.5 c/fm at t/t = 8. The transverse 
perturbation remains stable for the phenomenological parton model. 



VII. CLUSTERIZATION IN A MOLECULAR DYNAMICS CALCULATION 

When a perturbation of the zero order Vlasov equations is unstable, the linear approxi- 
mation breaks down for larger times. The development of the instabilities beyond the linear 
approximation can only be studied numerically. In the case of the mean field evolution 
(Vlasov equation), the numerical noise due to the finite number of test particles is a big 



problem [12]. It is very difficult to disentangle the numerical noise from the true fluctua- 
tions of the system acting as a seed for the development of the unstable modes. We have 
not attempted to solve this problem. Instead we study a different but closely related model. 
It is the molecular dynamical model with medium dependent masses. The local mean-field 
approximation to this model reduces to the Vlasov equation studied in the previous sec- 
tion. The model deals with classical partons whose masses are related to the density in the 
vicinity. The time dependent positions Xi(t) and momenta Pi(t) of the N partons follow 
Hamilton's equations 

dpi rrii 
—r — — v x mj 
dt Ei 

dxi L= Pi_ , s 

dt K J 

The parton mass is calculated from the gap equation with the scalar density in its vicinity. 

This requires the introduction of a range over which the density is calculated. We define the 

scalar density at position Xi(t) at time t by 

where the sum runs over all the particles j (with positions Xj(t), masses rrij(t) and energies 
Ej(t)) other than the parton i. It determines the mass of the parton % through the gap 
equation of the phenomenological parton model 
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^ = -*(*(«» ■ (42) 

The introduction of the finite range a (which we take a = 0.5 fm) has important conse- 
quences: (i) For a finite value of a, the mass of particle i at Xi and t cannot be determined 
by the masses of the particles j at the same time t, but only at earlier times (causality). As 
long as the masses vary slowly with time, the violation of causality may have no dramatic 
consequences, (ii) The size of the clusters into which the system breaks up will be strongly 
influenced by the value of a, but we have not yet studied this dependence. 

Using the gap equation for the phenomenological parton model, we have performed a 
simulation using 100 partons with initial temperature T = 180 MeV distributed in a fireball 



of radius r = 2 fm. The system expands and its density drops. Fig. [13] presents the evolution 
of the system at different times by showing the positions (x, y) of the partons projected on 
the z = plane. The partons tend to stay grouped in regions and structures of large local 
density. Eventually as the expansion proceeds the partons group themselves into compact 
clusters composed of at least two partons. There are no single partons. This is as close as 
one can do with the hadronization process in a model which does not know of color and not 
of quantum mechanics. The qualitative picture obtained from this calculation corresponds 
to the final stage of the dynamical break up of the system, which manifest itself first as a 
linear instability studied in the previous sections. The molecular dynamics evolution follows 
the system for times beyond the applicability of the linear response. 



VIII. SUMMARY AND CONCLUSIONS 

We have studied the instabilities of a plasma of partons. The plasma is treated in a 
quasi-particle model with an effective mass m, which depends on T and \i in the equilib- 
rium situation and on x, t for the expansion. The gap equation which relates the mass to 
the scalar density of the partons defines the underlying dynamics. We have studied two 
cases in parallel: (i) a gap equation derived from the NJL-Lagrangian in the mean field 
approximation and (ii) a phenomenological gap equation which is chosen so that certain 
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observables calculated from the lattice are reproduced by the quasi-particle model. While 
both models show a crossover transition for p — and a first order phase transition for 
p 7^ 0, only the phenomenological model has properties of confinement, in that the parton 
mass becomes very large for vanishing density. The study of these two, in some sense, com- 
plementary models, permits to get an idea on the model dependence of the results. The 
stability analysis has been performed on four levels: 

• For thermal equilibrium, one calculates the regions in the T — p B plane, in which 
the system is stable, metastable and unstable. These regions are calculated starting 
from the pressure as a function of T and p B . At thermal equilibrium instabilities 
are expected only at and due to the first order phase transitions. Since both models, 
the NJL and the phenomenological one show a first order phase transition for certain 
values of T and p B , the phase diagrams look rather similar from a qualitative point of 
view, but the values of T and p B at the phase boundaries are rather different. 

• Another interesting result found from the analysis of the phase diagram of the phe- 
nomenological parton model is the presence of the first order phase transition at finite 
baryon densities. Although we have followed the lattice QCD data at zero baryon 
density which show only a crossover transition arround T c the generalization to finite 
baryon density using a gap equation depending only on the scalar density leads to a 
first order phase transition. This simple model shows that the occurence of a first 
order phase transition in the highly excited matter with finite baryon density created 
in ultrarelativistic heavy ion collisions cannot be excluded from the lattice QCD data 
at zero baryon density. 

• The growth rates for the instabilities found in the thermodynamic analysis can be 
calculated from the Vlasov equation in a linear response analysis around thermal 
equilibrium. The neglect of the collisions in the linear response analysis in this work 
is justified when the instability growth rate is larger than the equilibration rate. We 
prove that for any quasi-particle model for which dm/dT < 0, one has a stable system 
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for pb = 0. In general the dispersion relation for frequency u and wave vector k of 
the instability has no scale and gives to oc k, since the gap equation corresponds to 
a zero range interaction. In the region of the mixed phase obtained from equilibrium 
thermodynamics, one finds solutions with Imu > 0, i.e., perturbations which grow 
exponentially in time. For the typical scale corresponding to the wave vectors k ~ 
1.0 fm _1 and as long as one is not too far away from the phase boundary^ the value of 
the growth rate Imu ~ 0.1 (fm/c) -1 is surprisingly small (compared to the expansion 
time scale of the plasma) . 

• The stability behavior far from equilibrium has been studied by investigating the 
Hubble- and Bjorken-solutions, (i.e., three- and one- dimensional expansion, respec- 
tively) of the Vlasov equation again by the method of linear response. The condition 
for the validity of this procedure is that the growth rate of the instability is faster 
than the expansion time (which itself must be larger than the relaxation time in order 
to be able to use the collisionless Vlasov equation). No instabilities are found for the 
NJL-model, while the phenomenological model quickly (r > 2r ) develops instabili- 
ties with significant growth rates (Imu ~ 0.5 to 1.0 (fm/c) -1 for k = 1.0 fm _1 ) even 
for ps = 0. We call these instabilities dynamical ones, since they are obviously not 
related to the first order phase transition but rather to the confinement properties of 
the phenomenological quasi-particle model. 

• The analysis of the instabilities based on linear response, be it at thermal equilibrium 
or for an expansion far from equilibrium, is limited to small times, if one finds expo- 
nential growth. The behavior over large times can only be studied numerically. We 
have attempted a first step in this direction, in that we have followed the expansion 
of a plasma described by the phenomenological quasi-particle model in a molecular 



1 The models based solely on the parton degrees of freedom may have not much physical meaning 
far in the confined phase. 
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dynamical calculation. We find indeed that the expansion of the system eventually 
leads to clusters with at least two partons. This is as far as one can get in a classical 
model which does not know color degrees of freedom. More work in this direction is 
intended. 
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FIG. 1. The temperature dependence of the parton mass in the two-flavor NJL model for [i = 0. 
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FIG. 2. The dependence of the parton mass on the chemical potential in the NJL model at 
T = 0. Note the expanded scale on the abscissa. 
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FIG. 3. The temperature dependence of the parton mass in the phenomenological parton model 
(PPM) at n = 0. 
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FIG. 4. The bag pressure density V as a function of the effective mass for the NJL model (dashed 
curve), and for the phenomenological parton model (solid curve). 
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FIG. 5. The dependence of the parton mass on the chemical potential in the phenomenological 
parton model at T = 100 MeV. 
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FIG. 6. The pressure density as a function of the baryon density in the NJL model for T = 0, 
20 and 45 MeV, solid, dashed and dash-dotted lines, respectively. The parts A-B and C-D of 
the T = curve correspond to the metastable states and the part B-C to the unstable state of 
the system. The straight line from A to D is the Maxwell construction for T = 0. A similar 
classification into unstable and metastable regions can be made for T = 20 MeV. For T = 45 MeV 
the system shows no first order transition. 
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FIG. 7. The solid line in the T — ps plane separates the stable region and the mixed phase. 
The mixed phase is devided into regions of instability ( "unstable" ) and metastability ( "ms" ) . The 
letters A to D at the bottom relate to T = and correspond to the points in Fig. |(| The results 
are for the NJL model. 
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FIG. 8. The pressure density as a function of the baryon density in the phenomenological parton 
model at different temperatures. The four curves (from bottom to top) correspond to temperatures 
T = 0, 50, 100, 120 MeV, respectively, while the insert magnifies the details around the origin for 
temperatures T = 100 and 120 MeV. 
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FIG. 9. The regions of the mixed, unstable and metastable phases in the T — pb plane in the 
phenomenological parton model. Note that the whole line ps = is part of the stable region. 
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FIG. 10. The contour plots for the instability growth rate F divided by the wave vector k in 
the T-pB plane for the phenomenological parton model. The solid line of this figure corresponds 
to the dashed one in the previous figure, which is the border of mechanical instability region. 
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FIG. 11. The parton mass as a function of the proper time for the three-dimensional homoge- 
neous expansion (solid line) and for the Bjorken expansion (dashed line). At t/tq = 1 one has a 
thermal distribution with T = 2T C . 
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FIG. 12. The instability growth rate T/k for a perturbation of the three-dimensional homoge- 
neous expansion (solid line) and for a longitudinal perturbation in the Bjorken expansion (dashed 



line) as a function of the proper time. Initial conditions like in Fig. 11 



33 



10 



- 



10 - 



E 



10 - 



- 



E 



■10 - 



t=0 fm/c 




10 



4 



i 1 1 1 r 

+ 



* 20 



+ 



10 10 

x(fm) 



10 10 

x(fm) 



FIG. 13. Expansion of a fireball of partons calculated within molecular dynamics. The positions 
of the partons with coordinates x and y are projected on the plane z = 0. Each cross represents a 
parton. At t = 20 fm/c there is no cluster with a single parton. 
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